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We present an analytical theory for the magnetic phase diagram for zigzag edge terminated hon- 
eycomb nanoribbons described by a Hubbard model with an interaction parameter U. We show that 
the edge magnetic moment varies as In U and uncover its dependence on the width W of the ribbon. 
The physics of this owes its origin to the sensory organ like response of the nanoribbons, demon- 
strating that considerations beyond the usual Stoner-Landau theory are necessary to understand the 
magnetism of these systems. A first order magnetic transition from an anti-parallel orientation of 
the moments on opposite edges to a parallel orientation occurs upon doping with holes or electrons. 
The critical doping for this transition is shown to depend inversely on the width of the ribbon. Using 
variational Monte- Carlo calculations, we show that magnetism is robust to fluctuations. Addition- 
ally, we show that the magnetic phase diagram is generic to zigzag edge terminated nanostructures 
such as nanodots. Furthermore, we perform first principles modeling to show how such magnetic 
transitions can be realized in substituted graphene nanoribbons. 
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Interest and activity in magnetic nanostructures has 
been driven by their possible application in nanoelec- 
tronic/spintronic devices, with graphene based systems 
grabbing a significant fraction of the attention. [1 The 
remarkable electronic properties of graphene [2] [3] with a 
Dirac like spectrum have made it a suitable candidate for 
many applications. [3H7] From a theoretical perspective, 
electron interactions and correlations effects on the hon- 
eycomb lattice contain many interesting phenomena[8 
including magnetism [9] and superconductivity. [TO] 

Magnetism at the zigzag terminated edges of graphene 
has been studied by first principles calculations [TTJ [12], 
and by a simplified effective Hubbard model [TOHTO] de- 
scribed by a hopping parameter t and a site-local repul- 
sion U. Ref. [TO] showed that the magnetism in graphene 
is robust to "shape disorder" of the nanostructure, while 
ref. [17] studied finite width graphene nanoribbons in- 
cluding the effects of doping. There are also studies of 
defect induced magnetism [20) |2T] and of magnetism of 
other nanostructures [T4] [22]. There are encouraging re- 
cent experimental signatures of magnetism[23-25 , along 
with suggestions [26] that extraneous effects such as re- 
construction would render the magnetism fragile. 

The origin of magnetic moment in zigzag edge ter- 
minated honeycomb nanostructures has been attributed 
to the edge states [27] [28] - localized electronic states 
which have most weight at the edges and die exponen- 
tially in the bulk. [29H5T] These states are of topolog- 
ical origin [32 and have been experimentally observed 
using scanning tunneling microscopy. [33 , 34 Magnetism 
at the edges is attributed to the Stoner mechanism(see, 
e.g., [35 ) and is best discussed in terms of a Landau 
theory. [35] The ground state energy of the system is ex- 



pressed as 

E(M) = a(U c - U)\M\ 2 + b\M\ 4 + . . . (1) 

where M is the magnetic order parameter, a, b > 
are positive constants that depend on the microscopies, 
and U c is a critical value of the on-site repulsion. For 
U < U c , the energy is minimum when M = 0, i. e., sys- 
tem is non-magnetic. At U = U c , there is a quantum 
phase transition to the magnetic state, and for U > U c 
one finds \M\ ~ y/U — U c . Stoner theory [36], based on 
linear response formulation, provides an expression for 
U c ~ g (l F ) •> wnere 9( £ f) is density of states of the bare 
system (U = 0) at the zero temperature chemical poten- 
tial sf- Thus for a zigzag edge terminated system, U c 
would vanish, since the density of states g(ep) diverges 
owing to the non dispersive nature of the edge states. 
The Stoner-Landau theory would therefore suggest that 
a zigzag edge will have spontaneous magnetization M for 
any U > 0, and furthermore that \M\ ~ yfjj for U <t. 

It is known [37 that a zigzag edge terminated nanorib- 
bon has a highly nonlinear response akin to that of sen- 
sory organs like eyes and ears. Their density response 
depends logarithmically on the magnitude of an edge po- 
tential applied at the zigzag edges. In this paper we show 
that such a Weber- Fechner response [38 of these nanorib- 
bons plays a central role in determining their magnetism. 
Our work leads to a "magnetic phase diagram" of lightly 
doped nanoribbons, including analytical expressions for 
the width W and U dependence of the magnetization, ex- 
citation gap, and the critical doping required to engender 
magnetic transitions. We also corroborate these results 
with variational quantum Monte Carlo calculations and 
show that the magnetism is robust to fluctuations. To 
the best of our knowledge, this is the first report that 
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FIG. 1. (color online) Two possible magnetic states of 
zigzag-edged nanoribbon, (a) antiferromagnetic (AF) and (b) 
ferromagnetic (FM). (c) Magnetic phase diagram of zigzag 
nanoribbons for U — 1.5 obtained from numerical mean field 
calculations. Inset: the critical doping required for AF to FM 
transition is inversely proportional to the ribbon width. 



clearly points out that the usual Stoner-Landau theory is 
inadequate to understand magnetism of hexagonal lattice 
nanoribbons. This gains added significance in view of re- 
cent developments of cold atom optical lattices, [39] where 
honeycomb lattices have been realized and studied. [40 

The simplest Hamiltonian that describes the energetics 
of interacting electrons in the honeycomb lattice is the 
Hubbard model 



H 



+ Ue ^2 n n n a - v^2 n icr (2) 



where c\ a is the operator that creates an electron of spin 
a at site z, n icr = c\ a Ci a is the number operator, t is the 
nearest neighbour hopping amplitude, U eo is the on-site 
Hubbard repulsion (eo = ^-t is a characteristic energy 
scale, U is dimensionless), and \i is the chemical poten- 
tial. The underlying triangular Bravais lattice has a lat- 
tice parameter a, and the width of the the zigzag edge 
terminated nanoribbon is denoted by W. 

Ground state of the zigzag edge terminated ribbon is 
obtained by a mean field analysis, where the four fermion 
interaction term is treated in the "magnetic channel" via 
the ansatz n^n^ \Ni ^ a n ia + MiSf - (\Nf - Mf) , 
where Sf = \ {n^ -n^), Ni = (^2 a n i<T ) is the mean 
occupancy, and Mi = (Sf) is the local magnetization 
at site i. The quantities Ni and Mi are to be deter- 
mined by enforcing the self consistency conditions of the 
mean field theory. Our detailed numerical calculation ex- 
ploits translational symmetry along the nanoribbon (see 
Fig. [TJa)) where we have used up to 4000 equally spaced 
points to sample the ID Brillouin zone. Fig. [T] shows the 
results for the values of physical parameters typical for 



graphene. We find that there are two possible magnetic 
configurations where the moments along the two edges 
of nanoribbon are oppositely aligned (the anti-ferro (AF) 
configuration, see Fig. [TJa)) and aligned in the same di- 
rection (the ferro (FM) configuration, see Fig.[TJb)). For 
both these configurations, the moments are concentrated 
on the sites at the edge layers - this point will be impor- 
tant in the discussion below. At half-filling (one electron 
per site), the ground state has an AF structure consis- 
tent with earlier results [11 . With the doping of holes 
denoted by S (S is defined as the doping per edge atom 
of the ribbon) , the ground state changes to the FM con- 
figuration at a critical doping [T7]. This critical doping 5 C 
required for the first order transition is dependent on the 
width of the ribbon and we find that S c ~ ^ as shown 
in Fig. [TJc), the "magnetic phase diagram". In the re- 
mainder of the paper, we develop an analytical theory of 
the physics behind this phase diagram. 

The continuum field theory [IT] that captures the 
physics of the Hamiltonian eqn. ([2| has the action (h = 1) 



s 



d 2+l 



Ue a 2 / d 2+ 



\%(r)(W T + U K )^a(r) 



r n(r) 2 -K(r)} 2 



(3) 



where r = (r,r), r = (x,y) is the position vector where 
^-coordinate is along the length of the ribbon and y along 
the width, r is the imaginary time that runs from to f3 
(inverse temperature), ^* = (ip* A +* ^B+a ^a-* ^b-*) 
is the array of Grassmann fields with A/B and +/— be- 
ing, respectively, sublattice and valley indices, n(r) = 
Eai/aCaW^aW is the number density, s z (r) = 
\ Ysava v^wairWaw M is the spin density, 1 is the 4 x 4 
identity matrix, and 
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with vp = eoa, and P — — iV, the momentum operator, 
r = T x e x + T y e y (t XjV - Pauli matrices in the sublattice 
space, G>x,y are spatial basis vectors). The action, written 
in a form that anticipates magnetism, can be studied by 
introducing a Hubbard- Stratanovich field m{r) to decou- 
ple the (s z ) 2 term in eqn. (|3|, while the n 2 term can be 
treated in a straightforward manner via a Hartree shift. 
The action becomes 

S[%m] = J d 2+1 rJ2K(r) (ld T + H K ) 



m(r)s z (r) + g.c. + TT ~ ^ 2 \m{r)\ 2 



1 



which upon integration of the fermion fields yields 

S[m) = InTTdet (— G~ x [m]) + — f d 2+1 r|m(r 

U e^a J 



(5) 
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FIG. 2. (color online) (a) Width dependence of magnetic mo- 
ment of the AF configuration for different values of U obtained 
from the full numerical solution compared with the analytical 
result eqn. |TT| ). Inset: The logarithmic dependence of the 
slope of Mq vs. a/W on U - points are obtained from the 
numerical simulations, solid line is from eqn. (11). 



Green's 



func- 



where the inverse 

tion is given by — G" 1 ^, r',m) = 
[lS r +HK + f{m(r)+m*(r)}l](5(r - r'). The 
magnetic ground state of the system can be described by 
a saddle point of the action eqn. (|6|. The saddle point 
field ra s satisfies the condition 
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We consider two different saddle point ansatzes for the 
undoped nanoribbons. For the AF configuration, we have 



mf F (y) = Ue a' 



M AF M AF 

S D (y-W) 6 D (y) 



(8) 

where M AF is the moment associated with one edge, and 
FM configuration has 
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where M FM is, again, the moment associated with one 
edge, and 5d is the Dirac delta function. Solution of 
M AF and M FM is aided by the observation that the 
right hand side of eqn. ^ is the magnetic moment den- 
sity response of a nanoribbon with applied edge Zeeman 
fields. As shown in ref. [37], this response is highly non- 
linear akin to that of sensory organs. Thus, for the AF 
configuration, eqn. ([7| reduces to 



M 



AF 



(10) 



For the FM configuration, M FM satisfies a similar equa- 
tion sans the factor of y/3 in front of U. For wide ribbons 



W > a, we find 
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(11) 
(12) 



in undoped case (hence the subscript 0). We further find 
that the AF configuration has a one-particle excitation 
gap given by 



e ~ W 1 6 i n 1 6 



(13) 



Fig. [2] presents a comparison of the analytical results for 
the AF configuration with the full numerical calculations. 
We find excellent quantitative agreement of the calcu- 
lated edge moment with the theory (eqn. ([TT])) over five 
decades of U (see inset [42] of Fig. |2|. Similar quantita- 
tive agreement is found for the energy gap eqn. (13). We 



also find excellent quantitative agreement between our 
theory and the numerical calculations of the FM state. 

Our theory predicts the AF configuration to be the 
ground state of the undoped ribbons for any width. This 
owes to the fact that the edge moment M AF is larger 
than M FM , 
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which, remarkably, is independent of U. 

The difference in the ground state energies 
(E AF — E FM ) per unit repeat distance along 
the length of the ribbon is estimated as ~ 

o e 

~ lnV3 a , ( (In y^) 2 , (6-7rUn a/3 , <T1\\ / a \ 2] 



-u 



-7r)ln>/3 i (U\\ ( g 

3 ttW 1 ^ 2 6 111 V 6 J J \ttW) 

where the leading term in a/W arises from the cor- 
relation energy (proportional to M 2 ), and the second 
term ~ (a/W) 2 also contains the kinetic energy of the 
electrons in the effective bands. The ground state of 
the undoped system is always the AF configuration the 
physics of which traces back to the sensory organ like 
response of zigzag-edge terminated ribbons. 

Upon doping the system with holes (S > 0) or electrons 
(S < 0), the edge magnetization changes. For the AF 
configuration, we find 



M (5) = M£ 



AF 



\6\ 



while, interestingly, for the FM configuration 

W. 



M FM (5) ~ Mq M 
where s(U) « f (l 



s(U)—S 2 
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FIG. 3. (color online) Dependence of C(U) on U. Critical 
doping 5 C for transition from AF to FM state is given by 
S C = C(U)§ (see eqn. p8l). 



to order a/W. We see that a first order transition from 
the AF to the FM configuration occurs at a critical dop- 
ing 



C(U) W , C(U) 



f ( 1+ RW6)) 



(18) 



Fig. [3] shows a comparison of this analytical result 
eqn. (18) with the numerical values of C(U), and again, 
excellent quantitative agreement is found over many 
decades of U. For larger values of U (such as that found 
in graphene), we find quantitative agreement up to about 
10%. This owes to the fact that larger values of U results 
in a small contribution from the bulk states, which is 
not captured in our edge mode based analytical theory. 
To the best of our knowledge, this is the first analyti- 
cal and quantitative theory of the magnetism in zigzag- 
edged ribbons. An important point brought about by 
the theory is that interpreting the magnetism in these 
systems via a Stoner-Landau theory of the form eqn. ([!]) 
may be too simplistic. Indeed, the energy functional has 
a non-analytic structure E(M) = AM 2 + BMlnDM 
(A,B,D constants), that has roots in the sensory organ 
like Weber-Fechner response of zigzag-edge ribbons. 

It is important to ensure that quantum fluctuations 
does not change the qualitative physics uncovered by the 
analytical theory. To this end, we performed variational 
Monte Carlo calculations of the ground state. Our trial 
wave function |^) = g D \FS,M), where D is the dou- 
ble occupancy operator, g is a the Gutzwiller factor that 
penalises double occupancy, \FS,M) is the filled Fermi 
sea state constructed by imposing an edge magnetization 
M on the edge layers. For the AF configuration M has 
opposite site on the two edges while for the FM config- 
uration M is equal on both edges. The optimal values 
of g and M are obtained so as to minimize the ground 
state energy. We have studied ribbons with W « 6a with 
U = 1.5. We find that the ground state at zero doping is 
the AF configuration, with an edge moment 0.20 which is 
expectedly smaller than the mean field value 0.27 owing 
to quantum fluctuations. Furthermore, the ground state 
changes to FM configuration at a critical doping S of 




FIG. 4. (color online) (a) Mean- field AF configuration in a 
zigzag edge terminated hexgonal nanodot (U = 1.5) (b) FM 
configuration obtained upon doping (5 = 0.1). 




FIG. 5. (color online) FM structure of the boron substituted 
zigzag edge terminated graphene nanoribbon obtained from 
first principles calculations. Dark atoms are borons and light 
ones are carbon. Hydrogen passivation of the edge atoms is 
also shown. 



0.025. Both the AF and FM edge moments are stable, in 
that they remain unchanged upon increasing the length 
of the ribbon (for a given width W), proving that the 
magnetism is robust. While these calculations are pro- 
hibitively expensive for the determination the full phase 
diagram, these results provide evidence for the correct- 
ness of our analytical theory. 

We have further investigated other zig-zag edge ter- 
minated nanostructures. For a undoped hexagonal nan- 
odot, we find that the magnetic configuration has an AF 
structure [14 as shown in Fig. Qa), which upon doping, 
changes to a FM type (see Fig. [4^b)). Doping in prac- 
tical applications can be achieved by gating. We have 
also explored possible chemical modification of graphene 
nanoribbons to produce an "internal doping" to engender 
the magnetic transition. Fig.|5]shows a FM structure of a 
zigzag edge terminated nanoribbon [43 with boron atoms 
substituted [44] in place of carbons, while the undoped 
nanoribbon has a AF configuration within the same cal- 
culation. As is evident, the results of this paper suggest 
many interesting possibilities of using zigzag edge termi- 
nated graphene nanostructures in applications. 

AM thanks CPDF, IISc for support. VBS is grateful to 
DST for Ramanunjan and MONAMI grants, and DAE- 
SRC for generous support. The authors thank Umesh 
Waghmare for help with the first principles calculations, 
R. Shankar and Jayantha Vyasanakere for discussions. 



* bsomnath@mrc.iisc.ernet.in 
^ amedhi@physics.iisc.ernet.in 



5 



* shenoy@physics.iisc.ernet.in 
[1] R. M. Westervelt, Science 320, 324 (2008). 
[2] A. H. C. Neto, F. Guinea N. M. R. Peres, K. S 
Novoselov, and A. K. Geim, Rev. Mod. Phys. 81, 109 

| (2009)| 

S. K. Pati, T. Enoki, and C. N. R. Rao, eds., Graphene 



[3] 
[4] 



[6] 



its fascinating attributes (World Scientific, 2011). 
J. Bunch, Y. Yaish, M. Brink, K. Bolotin, and 
P. McEuen, Nano Lett. 5, 287 (2005). 
[5] F. Miao, S. Wijeratne, Y. Zhang, U. C. Coskun, W. Bao, 
and C. N. Lau, Science 317, 1530 (2007). 
L. A. Ponomarenko, F. Schedin, M. I. Katsnelson, 
R. Yang, E. W. Hill, K. S. Novoselov, and A. K. Geim, 
Science 320, 356 (2008). 
[7] J. Wurm, M. Wimmer, I. Adagideli, K. Richter, and 
H. U Baranger, |New Journal of Physics 11, 095"022] 
| (2009)| 

[8] V. N. Kotov, B. Uchoa, V. M. Pereira, F. Guinea, 
and A. H. Castro Neto, ArX iv e-prints (2010), 
|arXiv: 1012.3484 [cond-mat.str-el] | 

[9] Z. Y. Meng, T C. Lang, S. Wessel, F . F. Assaad, and 
A. Muramatsu, Nature 464, 847 (2010)| 



[10] S. Pathak, V. B. Shenoy, and G. Baskaran, Phys. Rev. 
B 81, 085431 (2010) 



|11| Y. W. Son, M. L. Cohen, and S. G. Louie, |Phys. Rev. 
I Lett. 97, 216803 (2006)] 

[12J Y. W. Son, M. L. Cohen, and S. G. Louie, Nature (Lon- 
don) 444, 347 (2006). 
[131 N. M. R. Pere s, M. A. N. Araujo, and D. Bozi, |Phys.| 
| Rev. B 70, 195122 (2004)| 

[14] J. Fernandez- Rossier and J. J. Palacios, Phys. Rev. Lett. 
99, 177204 (2007). 



[15] S. Bhowmick a nd V. B. Shenoy, The Journal of Chemical 

Physics 128, 244717 (2008) | 

[16] O. V. Yazyev, Phys. Rev. Lett. 101, 037203 (2008) 



[17] 



J. Jung and A. H. MacDonald, Phys. Rev. B 79, 235433 
(2009). 



[18] T. Ma, F. Hu, Z. Huang, and H.-Q. Lin, Applied Physics 
Letters 97, 112504 (2010)] 



[19 



O. V. Yazyev, Reports on Progress in Physics 73, 056501 
(2010)| 

20J Y. Zhang, S. Tala patra, S. Kar, R. Vajtai, S. K. Nayak , 

and P. M. Ajayan, [Phys. Rev. Lett. 99, 107201 (2007)| 
[21] O. V. Yazyev and L. Helm, Phys. Rev. B 75, 125408 
(2007). 



[22] M. Ezawa, |Phys. Rev. B 76, 245 415 (2007) 
[23] H. S. S. Ramakrishna Matte, K. S. Subrahmanyam, and 
C. N. R. Rao, Ar Xiv e-prints (2009), |arXiv:0904.2739| 
[physics, chem-ph] 
|24| C. N. R. Rao, H. S. S. R. Matte, K. S Subrahmanyam, 

and U. Maitra, |Chem. Sci. 3, 45 (2012)1 
[25] V. L. J. Joly, M. Kiguchi, S.-J. Hao, K. Takai, T. Enoki, 
R. Sumii, K. Amemiya, H. Muramatsu, T. Hayashi, 
Y. A. Kim, M. Endo, J. Campos-Delgado, F. Lopez- 
Urfas, A. Botello-Mendez, H. Terrones, M. Terrones, and 
M. S. Dresselhaus, |Phys. Rev. B 81, 245428 (2010)| 



[26 

[27 

[28 

[29 
[30 
[31 

[32; 
[33; 

[34 
[35; 



J. Kunstmann, C. Ozdogan, A. Quandt, and H. Fehske, 
|Phys. Rev. B 83, 045414 (2011)| 

K. Nak ada, M. Fujita, G. D resselhaus, and M. S. Dres- 
selhaus, |Phys. Rev. B 54, 17954 (1996)| 
K. W akabayas hi, M. Fujita, H . Ajiki, and M. Sigrist, 

|Phys. Rev. B 59, 8271 (1999)1 

L. Brey and H. A. Fertig, |Phys. Rev. B 73, 1954 08 (2006) 
L. Brey and H. A. Fertig JPhys. Rev. B 73, 235411 (2 006) 
M. Zar ea and N. Sandler, |Phys. Rev. B 79, 16 5442 
(2009) 



S. Ryu and Y. Hatsugai, |Phys. Rev. Lett. 8 9, 077002 
(2002J] 

Y. Kobayashi, K. Fukui, T. Enoki, K. Kusakabe, and 

Y. Kaburagi, |Phys. Rev. B 71, 193406 (2005) | 

Y. Niimi, T. Matsui, H. Kambara, K. Tagami, 

M. Tsukada, and H. Fukuyama, Phys. Rev. B 73, 085421 

(2006). 

A. Altland and B. Simons, Condensed Matter Field The- 



ory (Cambridge University Press, 2006). 
K. Yosida, Theory of Magnetism (Springer- Verlag, 1996). 
S. Bhowmick and V. B. Shenoy, Phys. Rev. B 82, 155448 
(2010). 

See http:/ /en. wikipedia.org/wiki/Weber-Fechner Jaw. 
I. Bloch, J. Dalibard, and S. Nascimbene, |Nat. Phys. l^j 
267 (2012)| 

L. Tarruell, D. Greif, T. Uehlinger, G. Jotzu, and 

T. Esslinger, Natur e 483, 3 02 (2012). 

I. F. Herbut, [Ph7s. Rev. Lett 97, 146401 (2006)| 
The slope in the inset of Fig. [2] is the coefficient of the 
linear term obtained by fitting a polynomial M(a/W) = 
c + c 1 (a/W)+c 2 (a/W) 2 . 
[43] These calculations were done using Quantum Espresso 
code, [35] using plane- wave basis set and ultrasoft pseudo- 
potential. The energy cutoff for the plane- wave basis for 
wavefunctions is set to be 40 Ry. Electron exchange- 
correlation is treated with a local density approxima- 
tion(LDA, Perdew-Zunger functional). Nanoribbons are 
simulated using a supercell geometry, with a vacuum 
layer of 15 A between any two periodic images of the rib- 
bon. A k-point grid of 24x1x1 k points (periodic direction 
of the ribbon along x-axis) is used for sampling Brillouin 
zone integrations for this geometry. 
[44] L. S. Panchakarla, K. S. Subrahmanyam, S. K. Saha, 
A. Govindaraj, H. R. Krishnamurthy, U. V. Waghmare, 
and C. N. R. Rao, | Advanced Materials 21, 4726 (2009)| 
[45] P. Giannozzi, S. Baroni, N. Bonini, M. Calandra, R. Car, 
C. Cavazzoni, D. Ceresoli, G. L. Chiarotti, M. Cococ- 
cioni, I. Dabo, A. Dal Corso, S. de Gironcoli, S. Fabris, 
G. Fratesi, R. Gebauer, U. Gerstmann, C. Gougoussis, 
A. Kokalj, M. Lazzeri, L. Martin-Samos, N. Marzari, 
F. Mauri, R. Mazzarello, S. Paolini, A. Pasquarello, 
L. Paulatto, C. Sbraccia, S. Scandolo, G. Sclauzero, A. P. 
Seitsonen, A. Smogunov, P. Umari, and R. M. Wentzcov- 
itch, Journal of Physics: Condensed Matter 21, 395502 
(19pp) (2009)r 



